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Abstract 

We investigate relativistic spherically symmetric static perfect fluid models in the framework of the 
theory of dynamical systems. The field equations are recast into a regular dynamical system on a 3- 
dimensional compact state space, thereby avoiding the non- regularity problems associated with the Tolman- 
Oppenheimer-Volkoff equation. The global picture of the solution space thus obtained is used to derive 
qualitative features and to prove theorems about mass-radius properties. The perfect fluids we discuss 
are described by barotropic equations of state that are asymptotically polytropic at low pressures and, 
for certain applications, asymptotically linear at high pressures. We employ dimensionless variables that 
are asymptotically homology invariant in the low pressure regime, and thus we generalize standard work 
on Newtonian polytropes to a relativistic setting and to a much larger class of equations of state. Our 
dynamical systems framework is particularly suited for numerical computations, as illustrated by several 
numerical examples, e.g., the ideal neutron gas and examples that involve phase transitions. 
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1 Introduction 



The line element of a static spherically symmetric spacetime can be written as 

ds^ = -c^e'^^^Uf + e^^^^^dr^ + r^dO^ + ^ii^Od^^) , (1) 

where r is the radial area coordinate and c the speed of light, which we have chosen to introduce explicitly. 
A perfect fluid is described by the energy-momentum tensor T^jy = pu^u^ p^g^v + u^u^)^ where p is the 
mass-energy density and p the pressure as measured in a local rest frame; compatibility with static symmetry 
requires the unit vector to be parallel to the static Killing vector, i.e., = c~^e~^SQ. Einstein's field 
equations can be written so that the "relativistic potential" ^ is given by 



d^ _ 1 dp 

dr p -\- pc^ dr 



(2) 



while the remaining equations take the form 



dm ^2 /o \ 

— — = 47rr p , (3a) 
dr 

^l,.9I:lt(^ + JL)(^^^Iff)(^.?^y\ (3b, 

dr \ pC^ J \ mc^ / \ c^r / 

where G is the gravitational constant and where the mass function m(r) is defined by exp(2A) = (1 — ^^^^)~^ - 
To form a determined system, these equations must be supplemented by equations characterizing the properties 
of the perfect fluid matter, e.g., by a barotropic equation of state p{p). 

Eqs. ^ are the equations that govern hydrostatic equilibrium of a self-gravitating spherically symmetric perfect 
fluid body, Eq. (IHT^ is the famous Tolman-Oppenheimer-Volkoff (TOV) equation 19' T^. The TOV equation 
is fundamental to the study of relativistic stellar models, nonetheless it is associated with certain problems. 
One is that the equation is singular at r = 0, which leads to mathematical complications. It was only recently 
that such central issues as existence and uniqueness of regular static spherically symmetric perfect fluid models 
could be established Q]. A less known but perhaps even more problematic feature is that the right hand sides 
of Eqs. Q are not differentiable in p at p = for many barotropic equations of state, e.g., for asymptotically 
polytropic equations of state, discussed below. A similar problem occurs when p ^ oo; even though this 
limit is unphysical, this is still important. The relevance of the limit is due to the fact that there are several 
physical problems that are directly or indirectly governed by the limit p ^ oo. Consider, for example, the 
Buchdahl inequality 2GM/c^R < 8/9, where M is the total mass and R the surface radius, which implies that 
no star can have a larger gravitational redshift than 2 at the surface. All relativistic stellar models satisfy the 
Buchdahl inequality, however, the equality 2GM/c^R = 8/9, is associated with a solution of constant mass- 
energy density and infinite central pressure. Among other issues relying on differentiability of p{p) in the limit 
p ^ oo there is, e.g., the mass-radius relationships for large central pressures which we will discuss in the present 
paper. The derivation of the Buchdahl inequality exemplifies that sometimes certain rather unphysical solutions, 
associated with somewhat unphysical equations of state, are important for physical problems. Hence it is useful, 
indeed essential for certain problems, to understand the entire solution space associated with Eqs. including 
especially the limiting solutions with infinitely small or infinitely large pressures. In particular, in order to study 
generic behavior of solutions, one has to understand the entire solution spaces of large classes of equations of 
state. As indicated, the TOV equation is rather unsuitable for this task. 

The purpose of this paper is to introduce a formulation that avoids the above problems and allows us to (i) obtain 
a global picture of the entire solution space of a large set of equations of state and to (ii) probe relationships 
between the equation of state and global features like mass-radius properties. To this end we will develop a 
framework based on the theory of dynamical systems; Eqs. Q will be reformulated as a regular autonomous 
system of differential equations. A formulation that allows such a global treatment of solutions cannot exist for 
arbitrary barotropic equations of state. However, the large classes of equation of state that can be covered are 
of fundamental importance, as outlined in the following. 

Let us start by making some dimensional comparisons between General Relativity (GR) and Newtonian gravity. 
The dimensions relevant for self-gravitating perfect fluid models are mass, length, and time. In Newtonian 
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theory and in GR these are hnked through the gravitational constant G. Thus one can choose units so that 
one unit, e.g., mass, is expressed in terms of the other two. In relativity yet another fundamental constant 
appears - the speed of light c. This constant reduces the number of independent units to one, e.g., length. This 
restriction on the number of independent scales is of great significance since scaling laws are of great importance 
in gravitational theory, as in all branches of physics (see, e.g., |14j for references). 

The Newtonian limit of ^ and is obtained by defining the Newtonian potential^ and letting 

oo, thus setting the relativistic "correction" terms to zero. Note that p stands for the rest-mass density 
in the Newtonian case. Requiring that these equations be invariant under scalings of time and space forces 
a barotropic equation of state to be a polytrope p = where K > and n > is the so-called 

polytropic index. The appearance of c in relativity reduce the number of scale degrees of freedom from two 
to one since extra terms arise in the equations as compared to Newtonian theory; these terms originate from 
spacetime curvature and from pressure, which now also contributes source terms, and are proportional to c~^. 
Hence poly tropes do not leave the relativistic equations invariant under scale transformations. Only the linear 
and homogeneous equation of state p = Kp does. 

Invariance implies symmetry, and symmetry implies mathematical simplification. The existence of scale sym- 
metries allows one to choose all variables except one to be scale- invariant. Doing so leads to that the equation 
for the single non- scale- invariant variable decouples and thus one obtains a reduced set of equations for the 
scale-invariant variables. Hence poly tropes and linear homogeneous equations of state play a special role in 
Newtonian theory and relativity, respectively — they are the mathematically simplest perfect fiuid models 
these theories admit, and the associated equations can be formulated as symmetry-reduced problems.^ 

The central role these equations of state occupy is not due to gravity alone, but is also based on microscopic 
considerations of matter models. Again scaling laws can be used to motivate their special status, from a 
fundamental as well as phenomenological perspective. Consider, e.g., Chandrasekhar's equation of state for 
white dwarfs, where degenerate electrons contribute the pressure while the baryons supply the mass. In this 
case there are two limits — the non-relativistic electron limit, which leads to a polytropic equation of state with 
index n = 3/2, and the relativistic electron limit, which yields a polytrope with index n = 3. Another example 
is given by an ideal Fermi gas. For extreme relativistic fermions p = p(? describes the asymptotic behavior 
of the equation of state, where the fermions supply both the pressure and the mass-energy. 

This discussion suggests that it is natural to consider equations of state that are at least asymptotically scale- 
invariant, i.e., asymptotically polytropic or linear, and to classify equations of state according to their asymptotic 
limits; p ^ and p ^ oo. Their special status indicates that asymptotically polytropic/linear equations of 
state are the mathematically "simplest" equations of state — apart from their exact counterparts.^ 

In GR, asymptotically linear equations of state possess asymptotic invariance properties, however, surprisingly at 
first glance, also asymptotic poly tropes exhibit similar mathematical simplicity in the low pressure regime.^ The 
reason for this is the following: imposing spherical symmetry and requiring a static perfect fiuid equilibrium 
configuration results in a close connection between Newtonian gravity and GR. The spacetime symmetries 
prevent the existence of gravitational waves and gravitomagnetic effects, thereby reducing the number of degrees 
of freedom in GR to the same number as in the Newtonian case. Asymptotically, the relativistic equations 
of hydrostatic equilibrium coincide with their Newtonian counterparts provided that pj p(? when p 
0, and hence it follows that variables that are asymptotically scale- invariant in the Newtonian case are also 
asymptotically scale-invariant in GR.^ Accordingly, asymptotically polytropic equations of state can be naturally 
included in the relativistic formalism. 

The paper is structured as follows. In Sec. [21 the dynamical systems approach to relativistic stellar models is 
presented: by introducing asymptotically scale- invariant variables we reformulate the equations as a dynamical 
system. In Sec. |31 we discuss the assumptions on the equations of state that are needed for the dynamical 



-•^If one uses geometrized units c = 1, G = 1, there is no need for a distinction between <l> and 

^ Other equations of state may lead to greater simphfications for problems with extra symmetries, but it is only the polytropic and 
the linear equations of state that admit symmetry reductions in the general case in Newtonian gravity and relativity, respectively. 

^ Again we are referring to generic models, special cases, e.g., special static spherically symmetric models, might admit "hidden" 
symmetries, and may thus be mathematically even simpler or lead to explicit solutions, see |5|, il8i and references therein. 

^This is satisfactory since in this regime asymptotically polytropic equations of state are physically much more interesting. 

^This idea was first introduced by Nilsson and Uggla 27 , however, it will be exploited further here. 
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system to become regular. Thereby the notion of asymptotically polytropic/linear is made precise. The section 
concludes with a discussion about phase transitions. In Sec. 21 we investigate the state space of the dynamical 
system and study its global properties. The dynamical systems picture is subsequently translated into a physical 
picture in Sec. [SI In Sec. El we prove theorems concerning the relationship between the equation of state and 
global mass and radius features. To further illustrate our approach, we give several numerical examples in Sec. [3 
Finally, in Sec. 8 we give a concluding discussion, while we briefly compare the relativistic and Newtonian cases 
in the Appendix. 



2 The dynamical system 

Let us begin with some basic assumptions and definitions. Throughout this paper we will assume that the 
perfect fluid is characterized by a non-negative mass-energy density p and a non-negative pressure p, related by 
a continuous barotropic equation of state p{p) that is sufficiently smooth for p > 0.^ Define 

rN{p):='-'^, rip):=T-^\p) = P^, a{p):=A,, (4) 

p dp p dp pc^ 

where Fat is the standard Newtonian adiabatic index^ , however, it turns out to be more convenient to use 
the inverse index T, which in addition naturally incorporates incompressible fiuids (for which T = 0) into 
the formalism. The frequently used poly tropic index- function n{p) is defined via T(p) = n{p)/{n{p) + 1), or 
equivalently rAr(p) = 1 + l/n(p). 

Assumptions. For simplicity, although only necessary asymptotically in our framework, we assume that dp /dp > 
0. Then T > and p{p) is monotonic. If one impose the causality condition < c on the speed of sound 
Cs = y^dpjdp, then a < T. Below, further assumptions will be imposed by the dynamical systems formulation. 

A central idea of this paper is to obtain a regular dimensionless autonomous system from Eqs. ^ with a 
compact state space. To this end we first "elevate" r to a dependent variable and introduce ^ = Inr (r > 0) as 
a new independent variable. We then make a variable transformation from (m,p > 0,r > 0) to the following 
three dimensionless variables. 

The new pressure variable u;{p) is a dimensionless continuous function on [0, oo), strictly monotonically increasing 
and sufficiently smooth on (0, oo). We require that u;{0) = and a; ^ oo when p ^ oo, and hence a; > 0. The 
remaining freedom will be used to adapt co to features exhibited by the classes of equations of state of interest, 
cf. Sec. 01 The sign of u and v depends on whether m is positive or negative. In the following we restrict our 
attention to positive masses, i.e., we investigate a perfect ffuid solution only in that range of r where it possesses 
a positive mass-function.^ Moreover, we require that 1 — 2Gm/{rc^) > 0, i.e., we consider only solutions to Q 
that give rise to a static spacetime metric.^ Accordingly, we can assume u > and v > 0. 



^We will see in Sec. |21tliat, e.g., is sufficient, although this restriction can be weakened. Below we show how to handle even 
less restrictive situations like phase transitions. 

^However, note that p stands for the rest-mass density in the Newtonian case. 

^The region where m(r) < can be analyzed with the same dynamical systems methods that are going to be used in the 
following. The treatment turns even out to be considerably simpler, cf. the case of Newtonian perfect fluids 11 . 

^If a solution to (|3j satisfles 1 — 2Gm/{rc^) > initially at r = ro for initial data (mo,po), then this condition holds everywhere. 
This has been proved (for regular solutions), e.g., in 1 . Within the dynamical systems formulation this result can be established 
quite easily as we will see in Sec. ^ Solutions violating the condition 1 — 2Gm/ {rc^) > could be treated with the dynamical 
systems methods presented in this paper as well. However, we refrain from a discussion of such solutions here. 
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Starting from the transformation to the new variables u^v^uo yields the following system of equations: 

^=ii(3-ii + T/i) (6a) 

= y[-(l-T)h+{u-l){l + 2<jv)] (6b) 

^ = f^h , (6c) 
d In j9 

where h = h{u^ co) = — — = —v{l + cr)(l + au) , (6d) 

and / = /(oo) = dlnco/dlnp |p(a;)- Also, a and T are understood as functions of co. 
We now proceed by defining new bounded variables (/7, V, 1^) G (0, 1)^, 

l^U l^V l^UJ ^ ^ 

Introducing a new independent variable A according to dX/d^ = {1 — U)~^{1 — V)~^ yields the dynamical system 

^ = C/(l - [/)[(! - V){3 - 4U) -TH] (8a) 
^ = V{l-V)[{2U-l){l-V + 2aV) + {l-T)H] (8b) 

^ = -fn{i-n)H, (8c) 

where H = H{U, V, n) = {1 ^ a)V{l - U ^ aU) , (8d) 
and where T, <j, and / are now functions depending on Q. 

It turns out to be useful, indeed essential, to include the boundaries in our analysis so that our compactified 
state space consists of the unit cube [0, 1]^. To be able to discuss the dynamical system (|H| in a straight forward 
manner (e.g., to do a fixed point analysis) we require the system to be C^-differentiable on [0, 1]^. This is the 
case if Q{1 — Q)f{Q)^ ^{^)^ and are on [0, 1]. In the next section we will show that a broad class of 

equations of state satisfies these requirements. 



3 Equations of state 

In this section the consequences of the -differentiability requirement are examined in detail. There are three 
main building blocks: firstly, the requirement that 0(1 — Q)f{Q) is C^[0, 1] is shown to imply that the choice of 
uj{p) in is subject to certain restrictions; secondly, the assumptions on T are formulated in a precise way: 
T{ft) must be in an admissible pressure variable ft; this assumption defines the classes of equations of state 
we can treat in the dynamical systems framework; thirdly, a is made by the freedom in choosing cj. The 
section is concluded by illustrative examples. 

Part 1. The condition that Q{1 — Q)f{ft) is C^[0, 1] restricts the choice of the pressure variable uj{p) in 0. 
This condition is satisfied if and only if f{ft) is C^[0, 1] and C-^(0, 1), such that Vtdf/dVt and (1 — Q)df/dVt 
vanish in the limit ^ and ^ 1 respectively. Expressed in terms of /(p), a pressure variable uoij)) is 
admissible, if the conditions f{p) G C^[0,oo), limp^oo /(p) < oo, f{p) G C"'^(0,oo), and Imip^Q dhi f / dhip = 
limp^oo d\n f /d\np = are satisfied. It is not necessary, but convenient, to choose co in such a way that f{^) 
becomes strictly positive. Thus we require /o = /(O) > and /i = /(I) > in the following. -^^ 

-'^'-'Note that the Jacobian of the right hand side of contains d/dfl ^2(1 — ^2)], which when evaluated at ^2 = and Q = 1 

equals /o and /i, respectively. If these numbers are not zero the discussion of the dynamical system becomes easier. 
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Examples. The simplest example of an admissible pressure variable is provided by = {kpY ^ where a > is a 
constant and where /c is a positive constant carrying the dimension so that uj becomes dimensionless. In 

this example f(p) = dhiuo / dhip = which obviously satisfies the imposed conditions. Other simple examples 
include pressure variables described by uo{p) = {kp)^{±lnkp)^ for large p or small p, respectively. Here either 
a > and 6 G M, or a = and b <0. 

Part 2. Let us write the equation of state in the form 

kpc^ = ilj{kp) , (9) 

where /c > is a constant with dimension [p]~^; hence T{kp) = dlmlj{kp)/dln{kp). Recall that we assume 
dp/dp > 0, i.e., T > 0. 

Assumptions. We assume that there exists a variable Q (i.e., u;{p)) such that T{rt) is C^[0, 1]. Moreover, we 
require Tq = T(0) < 1 and Ti = T(l) = 1. 

Remark. The above assumptions restrict the equations of state to be asymptotically polytropic for low pressures 
and asymptotically linear at high pressures. Note that if one is only interested in classes of perfect fiuid solutions 
with bounded pressures, p < Pmax, then restrictions at = 1 are unnecessary, since then Q < l^max < 1- The 
case To = corresponds to an equation of state that is asymptotically incompressible when p ^ 0. Above we 
excluded the possibility Tq > 1, which includes the asymptotically linear case Tq = 1. This is not necessary, 
but models with Tq > 1 lead to solutions with infinite masses and radii, and are therefore not particularly 
interesting from an astrophysical point of view. Note, however, that T{p) might be > 1 for some range of p 
(corresponding to a negative polytropic index- function n(p)), which happens, e.g., for equations of state that 
cover the phenomenon of neutron drip. 

Let us now discuss the consequences of the assumptions for the equation of state. In a neighborhood ofQ = let 
us write T{Q) = To-\-Qv{Q)^ or equivalently, T(a;) = Tq-\-u;v{u;). Then the above assumptions can be expressed 
as follows: v{u;) is a continuous function, which is away from a; = 0, and satisfies lim^^o = 0-^^ By 

the definition of T we get dln{p~^^ p) / dlnp = uov^ which, using / = dhiuo / dhip^ gives rise to 



[ j r\^)v{io)dLo] =: sikpf" (1 + g{io)) , (10) 



kpc^ = s{kp)^° exp 

where 5 is a dimensionless constant. The function g{u;) is C^, and even for a; > 0, ^(0) = 0; moreover, 
uod^g/duo'^ vanishes as ^ 0. To first order we have g{uo) = f f~^{uo)v{uo)dix). 

Analogously, in a neighborhood of 1] = 1, one can show that the assumptions can be translated to 

kpc^ = (71 - l)-\kv) {l+~g{Lo-^)) , (11) 

where 71 > 1 is a dimensionless constant, satisfying 71 < 2 if we require (asymptotic) causality. Hence we 
have shown that the above assumptions on T are satisfied if and only if there exists a pressure variable co and 
functions ^, g with the above properties, such that the equation of state Q can be written in the asymptotic 
form and (|TT|) . respectively. 

Examples. For a simple example consider the pressure variable uj = {kp)^ and take g{uj) = const uj. Then (PTTj) 
reads kpc^ = {kp)^^ (1 + const{kpY). 

Remark. If the assumptions hold, then it is always possible to choose uo in such a way that dT/dft vanishes as 
1] ^ and 1] ^ 1. We can even achieve 

g=0(f]'=) {n^O) and ^=0{{l-nf) (fi^l) (12) 

for an arbitrary k > 0. This is based on the fact that if co is an admissible pressure variable then so is a; = cj^ 
(0 < / < 1)^2 ^ dT/dn = O(r2(i-^)/0 for n^O, and analogously for Q ^ 1. 



■"^ -"^ Compare also with the discussion about admissible (jo(p). 
^^Its logarithmic derivative / = d\nuj/dp satisfies / = 
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Part 3. The third imposed condition is that a{Q) is C^[0, 1]. Firstly, the assumptions on T clearly imply that 
a{Q) is continuous with ctq = cr(0) = and ai = a{l) = 71 — 1. Secondly, a{rt) is C-^[0,1], if ft is chosen 
appropriately: a straight forward computation of dcr/dVt yields da / dVt = (j[Vt{l — Vt)]~^ {1 — T) . For 1] — > 1 
we obtain lim^^^i dcr/dQ = — (71 — l)/;f ^ lim^^^i dT /d{l — Q)^ hence lim^^^i dcr/dfl exists by the assumptions on 
T. In the limit 1^ ^ we obtain limj^^o diner /d\nQ = \imQ^o{Q/cr) dcr/dQ = /^^{I — Tq). If /q""^(1 — To) > 1, 
then limj^^o da/dQ exists and lim^^^o dcr/dVt = 0. If /q~^(1 — Tq) < 1, then da/dft ^ 00 for ^7 ^ 0.^^ Therefore, 
provided that Q is chosen so that /o < (1 — Tq) = 1/(1 + ^o), we conclude that G C^[0, 1]. In addition, in 
analogy with the discussion involving T, we can always achieve 

^ = 0{n'^) (n^o) and ^ = o{{i-n)'^) (n^i) (13) 

for arbitrarily large k. 

Remark. Consider the equation of state Q. We choose cc; as a function of the dimensionless kp as already 
indicated above, i.e., u — uj{kp). Then the equation of state can be written in the implicit form 

kp = (l){uj) , kpc^ = x{^) , (14) 

where (j) and x contain an arbitrary number of dimensionless parameters. The derived quantities /, T, 
and a then read f{uj) = {dln(j)/dlnuj)~^ ^ ^{^) = {dlnx/dlnuj) {dln(l)/dlnuj)~^ , and cr(cj) = (j){uj) / x{^) - Note 
that /(cc;), T(cc;), a{uj) are functions independent of the dimensional parameter k; therefore the whole class of 
equations of state (|T^ . or equivalently ©, parameterized by k is described by one dynamical system with 
specified /(l^), T(l]), and a{Q). 

Example. Relativistic polytropes. The relativistic generalization of the Newtonian adiabatic index Fat (as 
defined in (^) is the adiabatic index 



Clearly, F reduces to F^r when c 00.^^ Regarded as a function of the adiabatic index F becomes a 
C^[0, 1] function for the considered class of asymptotically polytropic/asymptotically linear equations of state. 
In particular, F(0) = FAr(O) = Tq "'^ = 1 + 1/no, and F(l) = 1 + <ji = 71. Whereas F^r = const defines the 
polytropic equations of state, F = const gives rise to the so-called relativistic polytropes, usually written as 

pc' = {p/K)'^^c^+p/{r-l), (16) 

where i^T is a positive parameter. Requiring dp /dp > and causality of the speed of sound Cs = ^^/dp/dp < c 
implies 1 < F < 2, where F = 2 corresponds to the case Cg ^ c when p ^ 00. The class ((16|) of equation of 

state can be expressed nicely in the form Q, if we take k = (c^(F — l)i^~^/^) : 

kpc^ = (j^) ((M'/^ + kp) . (17) 

Obviously, the relativistic polytropes are of the form (|Tnj) and (|TT|) by the simple choice uj = {kpY with a < 
(F — 1)/F, whereby T(yt) and cr(yt) are C^[0, 1]. To simplify the expressions for T(yt) and cr(yt) as much as possible 
we take uj = {kp)^^-^^^, i.e., / = (F - 1)/F, which results in a{n) = (F - l)n and T{n) = (1 + (F - l)n)/r. 
Hence the right hand side of the dynamical system consists of pure polynomials. Note that only one 
(dimensionless) parameter appears in the dynamical system, namely F, i.e., the entire one-parameter family 
parameterized by k is represented by a single dynamical system, as expected. 

Example. Ideal neutron gas. The equation of state of a degenerate ideal neutron gas is implicitly given by 

^P=^ (^(1 + - 1) + + (1 + ^')'^']) (18a) 

kpc^ = ^ (x(l + x^y^\l + 2x^) - \n[x + (1 + x^)!/^]) , (18b) 



-•^^If fo = {l — To), then in many cases liniQ^o da/dfl exists. Note, however, that this does not hold in general. 
■'^^Recall, however, that in the Newtonian case p stands for the rest-mass density. 
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where x is essentially the Fermi momentum. A straight forward way to treat (|T^ is to choose co = 
whereby Q = x/{l -\- x): Firstly, x is an admissible pressure variable, as f{^) is a sufficiently smooth function 
satisfying / = 1/5 + 0{n^) {n 0) and / = 1/4 + 0((1 - Q)^) 1). Secondly, a{n) is smooth with 

a = Vt^/b + 0{Vt^) and a = 1/3 + 0((1 - Q)'^) as 1] ^ 0, 1 respectively. In particular, cji = 1/3, i.e., 71 = 4/3. 
Thirdly, T{Q) interpolates smoothly between To = 3/5, i.e., no = 3/2, and Ti = 1. Moreover, dT/dQ = 0{Q) 
and dT/dQ = 0(1 — Q) for ^7^0,1. Therefore, the ideal neutron gas is easily described in the dynamical 
systems formalism, moreover, our formalism also covers more sophisticated equations of state that are used to 
model neutron stars. 

To conclude this section we outline how phase transitions and composite equations of state can be incorporated 
into the dynamical systems framework. Consider an equation of state p = p{p) that is piecewise continuous on 
[0, 00) and piecewise sufficiently smooth on (0, 00). This implies that at certain values pj of the pressure, p{p) 
makes a jump or a kink, corresponding to a phase transition of ffist or second kind; T(p) and a{p) exhibit similar 
behavior. The associated "weak solutions" p{r) and m(r) of equations © are required to be continuous, while 
p{r) is continuous only if p{p) is. Also the variables v and /7, V are continuous functions if p is; if p jumps 
from some value pi to p2 at pj, then u^v make a jump according to Ui/u2 = Vi/v2 = pi/ P2-, and analogously 
for U^V. Choosing a conventional continuous pressure variable we obtain piecewise continuous/smooth 
orbits, monotonic in 1], that have possible jumps in /7, V at Vtj = ^(jpj). 

For practical purposes a slightly different viewpoint turns out to be more convenient. Assume for simplicity 
that there is only one jump of at pj. Choose two smooth equations of state pi{p), Pii{p)i such that p = pi 
for p < Pj and p = pu for p > Pj. Then the equation of state p is simply obtained by switching from pi to pu 
at Pj. Correspondingly, we deal with two different state spaces, in particular we can have two different pressure 
variables Qj and Qu. /^From this point of view, a fluid solutions associated with p starts in one state space, 
and ends in the other state space; the jump in between appears as a map from one state space to the other; 
besides the jump in [/, V described above it entails a jump in 1], from ^i{pj) to ^ii{pj). An explicit example 
will be given in Sec. 



4 Dynamical systems analysis 

In this section we study the dynamical system ((Hj); in particular we investigate the global dynamics. The arena 
for the analysis is the state space, i.e., the unit cube [0, 1]^, which is endowed with different vector fields and 
associated flows that depend on the equation of state. 

Central to the dynamical systems analysis are the equilibrium points of the system. The fixed points as well as 
their associated eigenvalues^^ are listed in Tabled] Note that the equation of state enters only via its asymptotic 
properties: both the location of the fixed points and the eigenvalues depend only on no and 71. The state space 
and the fixed points are depicted in Fig. 

We observe that all six faces of the cube are invariant subspaces. On each of the four side faces, the induced 
dynamical system possesses a simple structure; the orbits on the side faces^^ are depicted in Fig. [2l 

The dynamical system on the 1] = subset is in fact Newtonian, cf. Appendix IXI since a = on 1] = 0, 
and V coincide with Newtonian homology variables and thus the system represents the Newtonian homology 
invariant equations for an exact polytrope with index no- On the = subset the flow of the dynamical 
system exhibits rather complicated features. In particular, no serves as a bifurcation parameter: there exist 
bifurcations for no = 0, no = 3, and no = 5, as indicated by the properties of the fixed points. The case no = 5 

-•^^Namely, x = pF/mnC, where mn is the neutron mass. The constant k is given by k = X^/mnC^, where An = h/mnC, for details 
see, e.g., |^, p. 23ff. 

■•^^The eigenvalues of the linearization of the system at the fixed point (together with the eigenvectors) characterize the fiow of 
the dynamical system in a neighborhood of the fixed point (Hartman-Grobman theorem). For an introduction to the theory of 
dynamical systems, see e.g. 9 . 

■•^^Note that whereas the orbits in the interior of the cube represent perfect fiuid solutions, the orbits on the side faces cannot 
be interpreted in physical terms. This is essentially because the variable transformation ^ cannot be inverted for, e.g., u = or 
V = 0. In a modified dynamical systems formulation, where lj is replaced by a variable depending on <l>, certain side faces may be 
interpreted as representing vacuum solutions. 
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Fixed point 


U 


V 


n 


Eigenvalues 


Restrictions 




1 







1,1,0 




Lo 


3 
4 







_3 1 
4 ' 2 ' 




u 








^0 


3,-1,0 




u 


Uo 


1 





0, -{l-Uo) , -/o(l-^o) 


no = 


Bi 





no + 1 





no-3 1 /o(no + l) 




no+2 


no+2 ' no+2 ' no+2 




B2 





1 





no 1 f 
no + 1 ' no + 1 ' -^^ 




Bs 


1 


1 





0,0,0 




B^ 


no— 3 


2no+2 





^{5-no±iVb) , -/o/3(no^-l) 


no > 3 


2(no-2) 


3no + l 


Ti 





1 


1 


-71 , 2(71 - 1) , /i7i 


1< 71 < 2 


T2 


1 


1 


1 


71(71 - 1) , -2(71 - 1) , /i7i(7i - 1) 
-:^t(37i-2±Zv^) , :^/i7i' 


1< 71 < 2 


Ts 


1 

2 


2 


1 


1< 71 < 2 



Table 1: Local properties of the fixed points. We use the abbreviations b = —1 — 22no + 7no, /3 = 
and c = -71^ + 4471 - 36. 




Figure 1: The compact state space for the dynamical system with equilibrium points. The position of Bi 
and ^4 depends on no (^4 is only present for no > 3) and the position of T3 depends on 71. 



is exactly solvable: the orbits can be represented implicitly by 

$ = VUV^^^{U{S - 9 V) + 7 V - 6) - d(l - Uf^^il - Vfl^ = , (19) 

where the parameter d is required to satisfy d > —9/8. The fixed point Bs is represented by d = —9/8; 
— 9/8 < (i < defines a one-parameter family of closed orbits, which we will henceforth denote by C; d = 
characterizes the orbit that connects L2 and Bi; the remaining orbits correspond to positive values of d, 
cf. Fig. |2(i)| Typical orbits for a selection of values of no are depicted in Fig. [3 cf. also the analogous non- 
compact figures in the standard literature, e.g., p. 201 in 20 . 

The dynamical system on the Q = 1 subset can be regarded as describing the relativistic scale-invariant equations 
for linear equations of state p = (71 — l)pc^, a fact that refiects itself in the decoupling of the ([/, V)-equations 
from the system for purely linear equations of state. Typical orbits are depicted in Fig.I21(cf. also, e.g., J^, 
where an analogous non-compact phase portrait is given). 

The global dynamics of the dynamical system is determined by the fixed points on the boundaries of the 
state space, which enables us to describe the asymptotics of the interior orbits. The statement is made precise 
in the following theorem. 
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(a) U = 



{h)U = l 




{d)V = l 




(g) 3 < no < (11 + (h) (11 + 8v^)/7 < 

8V2)/7 no < 5 




(j) no > 5 



(k) 71 = 1.1 



(1) 71 = 1-9 



Figure 2: In (a), (b), (c), and (d) orbits for the four side faces are shown. The U = 0, U = 1, and V = 1 faces 
involve a and T, which depend on the equation of state. For the plots the relativistic polytrope with F = 4/3 
has been chosen; however, the qualitative structure is the same for all equations of state. In (e), (f), (g), (h), 
(i), and (j) orbits on the Newtonian subset Cl = for different no are presented. The subfigures (k) and (1) show 
orbits on the Q = 1 subset for some values of 71. 
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Theorem 4.1. The a-limit^^ of an interior orbit is a fixed point on Li, L2, or the fixed point T3. The uo-limit 
of an interior orbit is always located onVt = and it is one of the fixed points Bi, B2, or when no 7^ 5. // 
no = 5; then the uo-limit can also be an element of the 1-parameter set of closed orbits C . 

Proof. The proof of the theorem is based on the monotonicity principle ^Sl • if there exists a function Z (defined 
on a closed, bounded, future (past) invariant set S) that is monotonic along the fiow of the dynamical system 
X = /(x), then uj{S) {a{S)) is contained in {x G S\Z = 0}. 

The function Q is such a monotone function on the state space, i.e., dQ/dX < 0, except on 1] = 0, 1^ = 1, and 
V = 0. It follows that the cj-limit of an interior orbit is located on 17 = or y = 0, while the a-limit of an 
interior orbit must lie on 17 = 1 or 1/ = 0. However, it follows from the known orbit structure on V = and 
the local fixed point analysis that V = cannot possess an cj-limit. Similarly, it follows that only Li and L2 
constitute a-limits for interior orbits. 

It remains to determine the cj-limit (a-limit) sets for interior orbits on 17 = 1 and 17 = 0. Again we make use 
of the monotonicity principle, however, to simplify the discussion we use the uncompactified variables {u^v). 

On (7 = 1, a monotone function is given by^^ 

Z3 = [-(71 - l)hi{u,v) + 371 - 2]Huv)-P'{l + 2(71 - l)v)P'-^ , (20) 
where the non-negative function hi and the non-negative constant pi are defined by 

4(71 - 1)2 

hiiu, v) := 71(1 + (71 - l)u)v = -h\n=i , Pi = ^ ^ . (21) 

71(571 - 4j 

Its derivative along orbits reads 

^ = - Y (371 - 2)(2 - /ii)2[(7i - l)h^ + 371 - 2]-^Zs . (22) 

Application of the monotonicity principle yields that only T3 is a possible a-limit set on ^7 = 1. 

We now turn our attention to the c<;-limit sets of interior orbits, which we have shown to lie on 17 = 0. For 
no < 3 the function Zi = uv^ is monotone along orbits, dZi/d(_ = Zi{2u + '^(3 — 4To)). For no > 3 define 



{uv) 



2/(no-l) 



no — 5 2v v'^ uv 



{no - 1)2 (no + l)(no - 1) 2(no + 1)^ (no + 1)^ 



(23) 



and observe that dZ2/d^ = (no — 5) (no + l)~'^{no — l)^((ii'u)2/(^o~^) (2(no + 1) — (no — 1)'^)^, i.e., Z2 is monotonic 
on orbits for 3 < no 7^ 5. (Note that the fiow is non- vanishing on v = 2 (no + l)/{no — 1) except for at B4.) 
Thus, by the monotonicity principle, taking also the solution structure on the boundaries into account, only 
the fixed points 5i, B4 can serve as C(;-limit sets on 17 = 0. The exceptional case is no = 5: then there also 
exists a 1-parameter family, C, of closed orbits (periodic solutions) described by (see also (|T^ ). which acts 
as an cj-limit set for interior solutions. ■ 

Theorem 4.2. //O < no < 3 (no = 0)^ then all orbits end at B2 (^4). 



Proof We distinguish two cases, no < 3 and no = 3. For no < 3 the proof is trivial; each interior orbit converges 
to B2 for < no < 3 (to L4 when no = 0), which follows from the proof of Theorem l4.1l and the local dynamical 
systems analysis (Table QJ. 

In the case no = 3, is not hyperbolic. Nevertheless, by applying center manifold theory, we show that no 
interior orbit converges to Bi as A ^ 00: 

It is advantageous to investigate the problem in adapted uncompactified variables {x > 0, z > 0}, defined by 

u = x{no - 2)/{no ^ 1) v = {uq ^ 1) - x ^ y uo = z , (24) 

^^For a dynamical system x = /(x), the a;-limit set (a-limit) of a point x is defined as the set of all accumulation points of the 
future (past) orbit of x. Limit sets thus characterize the asymptotic behavior of the dynamical system. 
^^This monotone function was found through Hamiltonian methods, see |^ and chapter 10 in 13 . 
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The fixed point Bi is represented by (0, 0, 0) and the dynamical system © reads 

x= + Nx{x,y,z) 

y= y + Ny{x,y,z) (25) 

z= -fo{l^no)z + N^{x,y,z), 

where N^^ Ny, denote the nonlinear terms. is given by 

= x{3 - 4T{z) - x/4 + (x - y)T{z)) + Ci{x, y, z)a{z) + C2{x, y, z)a\z) . (26) 

A center manifold of the system (|25j) . i.e., an invariant manifold tangential to the center subspace at the fixed 
point, is represented by {(x, hy{x)^ hz{x))\0 < x < e} satisfying 

dxhy{x) x{x, hy{x), hz{x)) = y{x, hy{x), hz{x)) dxh^ix) x{x, hy{x), hz{x)) = i(x, hy{x), hz{x)) , (27) 

and the tangency conditions /i(0) = 0, dxh{0) = 0. We obtain 

hy{x) = -x^/2-x^ ^0{x^) hz{x)=0. (28) 

Note in particular that the center manifold lies in the Newtonian subset of the state space. The center manifold 
reduction theorem states that the fiow of the full nonlinear system is locally equivalent to the fiow of the reduced 
system (see, e.g., |2j) 

X = Nx{x, hy{x), hz{x)) = x^ 12 + 0{x?) 

y = y (29) 
z = -/o(l + no)^. 

The fiow given by ((2^ clearly prevents interior solutions {x > O^y^ z > 0) from converging to 5i, which proves 
the claim. ■ 

Collecting the results of the local and global dynamical systems analysis yields the following. Interior orbits 
originate from Li (a two-parameter set); L2 (a one-parameter set); and T3 (a single orbit). They end at Bi 
(when no > 3; a one-parameter set); B2 (a two-parameter set); I/4 (when no = 0; a two-parameter set); ^4 
(when no > 3; a single orbit when no < 5 and a two-parameter set when no > 5); and at the periodic curves C 
(when no = 5; a two-parameter set). 

An orbit in the state space corresponds to a perfect fiuid solution. Note that the perfect fiuid solution is only 
represented in that range of r where p > 0, i.e., the fiuid body, and where m > 0; cf. Sec.[2l If the configuration 
has a finite radius then the fiuid solution is joined to a vacuum solution, a Schwarzschild solution, at r = i?. 

The one-parameter set of re^i^/ar perfect fiuid solutions, i.e., solutions of with a regular center of spherical 
symmetry, appears in the state space as the one-parameter set of orbits that originate from L2. When dp /dp > 
is assumed, p{r) < p{r) holds for regular solutions, where p = 3m/(47rr^) is the "average density". Translated 
to the state space variables we obtain U < 3/4. In order to show this inequality within the dynamical systems 
framework, consider the surface U = 3/4 and the fiow through this surface: dU/dX\ij=s/4 = — 3T(1 + + 
3<j)y/64 < 0, since T > when dp/ dp > 0. Hence, /7 = 3/4 acts as a "semi-permeable membrane" for the fiow 
of the dynamical system. By noting that the unstable subspaces at I/4 are characterized by U < 3/4 the claim 
is established.^^ 

If a perfect fiuid solution satisfies 1 — 2Gm/{rc^) > initially at r = ro for initial data (mo,po), then this 
condition holds everywhere. For a proof (involving regular solutions) see, e.g., |1|. In the state space picture 
this result can be established quite easily. By construction, the initial data is a point in the interior of the state 
space and the orbit passing through this point represents the associated perfect fiuid solution. Since V < 1 for 
the orbit, we obtain that 

2Gm 1-V 

^^We observe that in the incompressible case, p = const (i.e., T = 0), the regular solutions are located on U = 3/4. The equations 
are explicitly solvable in this case. Note, however, that we have to cut the state space at some value of 17, because T = is not 
compatible with the requirement Ti = 1 to obtain a. dynamical system up to and including fl = 1. 
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holds everywhere within the fluid body. We can even show that 1 — 2Gm/ (rc^) > const > for all r: the 
considered orbit can have three possible end points; suppose that the orbit ends in the fixed point which is 
the only non-trivial case since then V ^ 1 when A ^ oo. We will see in the subsequent section that such an 
orbit gives rise to a perfect fluid solution possessing a flnite radius R and a total mass M, which are related by 
GM / [Rc^) = B/(2B + 5), where B and s are positive numbers. Hence the claim is established. 

By introducing the function 

-^~2{2M + V) ' (^^^ 

where 

we can write the relativistic Buchdahl inequality as 

^1 > . (33) 

Note that the surface inequality < 4/9 is obtained by setting 7^ = 0. The flow on the surface ^1 = is given 
by 

d^i _ G{^-AU){l-U^(jUfV 

li\ ~ (2(1-/7)+ (7/7)2(3(1 -/7)+a/7) ^ ^ 

and hence, since U < 3/4 for regular solutions, it follows that the derivative is positive (zero) when U < 3/4 
{U = 3/4). The Buchdahl surface ^1 = and the U = 3/4 surface are depicted in Fig. 01 for the relativistic 
polytropes with F = 3/4 and F = 2 (note that the Buchdahl-surface is affected by the equation of state through 
a{Q)). Since regular solutions start at ^ = and since the flow on the Buchdahl surface is directed toward 
y = 0, the regular solutions cannot pass through the Buchdahl-surface. Hence, for regular solutions only a 
limited part of the state space is accessible (cf. Fig. This constitutes the dynamical systems proof of the 
Buchdahl inequality. 




5 Translating the state space picture to a physical picture 

In this section we translate the state space picture to the common physical variables. To this end recall, flrstly, 
that an orbit in the interior of the state space stands for a fluid solution {m{r) ^ p{r) ^ p{r)) (r > 0) inside the 
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fluid body via the transformation Secondly, consider the following auxiliary equations: 



dr 



(l-C/)(l-y)r 



{l + a)V{l-U + aU)p , 



dX 
dp 

d\ 

Gm 

~ l-V + 2aV 



dm 



C/(l - V)m 



da 
dX 



-(1 + cr)(l - T)V{1 -U + aU)a 



aV 



U 



AttG \1-U J \1-V + 2aV 



aV 



m 



U 



47rG3 \1-U J \1-V^ 2aV 



aV 



(35a) 
(35b) 
(35c) 

(35d) 

(35e) 



Based on an equation of state in explicit form 
replace in the last two equations by 



kc^ 



, i.e., kpc^ = ip{kp)^ and a pressure variable uj = u{kp) we can 



p-^ = kc^a{n) 



{kp){Q) ' 



(36) 



where {kp){n) = (^^/^/(l - ll)^^ in the case u = {kpY . 

Equations (|35a|) can be used to obtain an intuitive picture of where solutions gain mass and radius in the state 
space. Since these equations are independent of 1^, it is possible to visualize dhim/d\ and dhir/d\ as contour 
plots on the /7, V plane, see Fig. ^ The picture suggests that no mass is acquired near the side faces /7 = 
and V = 1, and hence the fixed points Bi and B2 (and I/5, if present) are attractors for solutions with finite 
masses. In contrast, dhir/d\ does not vanish near 5i, but only near B2 (and I/5, if present). Therefore, B2 
is the only attractor for solutions with finite radii. Solutions converging to B^, (or to a periodic orbit C when 
no = 5) for A ^ 00 acquire both infinite masses and radii. 




Figure 4: Contour plots of and The white regions correspond to the minimum value of the functions 
(= 0) and the black to the maximum value (=1). 

To obtain a rigorous derivation of the asymptotic mass and radius features of the various solutions we must 
combine equations ((35j) with the local dynamical systems analysis in the neighborhood of the fixed points (and 
the periodic orbits C when no = 5). We only state the results here (for a detailed discussion of the corresponding 
and rather similar Newtonian situation, see PTI): 

Solutions that originate from: (i) Li are solutions with negative masses for small r; when they have acquired 
sufficient positive mass, so that the mass function m has become zero, they appear in the state space through Li, 
cf. the Newtonian case discussed in pjj); (ii) L2 have a regular center and are therefore solutions of particular 
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physical interest; (iii) T3 have a non-regular center (p ex for r ^ 0) and describe the limit when the 
central pressure goes to infinity. One can show that T3 is associated with a self-similar solution admitting a 
spacetime transitive symmetry group. It is a special case of a solution found by Tolman JL^, but has been 
rediscovered many times, see and references therein. 

Solutions that end at: (i) Bi (no > 3) have infinite radii but finite masses; (ii) B2 and (iii) L4 (no = 0) have 
finite masses and radii and are therefore the most interesting solutions; (iv) B4 (no > 3) and (v) C (no = 5) 
have infinite radii and masses {m{r) ~ /^(no-3)/(no-i) ^^j^^ ^.^^^ 

Let us now take a closer look at the 2-parameter family of orbits that converges to ^2 as A ^ 00. All such 
orbits correspond to perfect fluid solutions with finite radii R and masses M. To describe the behavior of the 
physical observables as r ^ i?, we will use the local dynamical systems results together with equations Q, 
13), or Recall that the solutions of the dynamical system with specified T{Q) and g(Q) represent the 

totality of perfect fiuid solutions corresponding to a one-parameter family of equations of state, cf. ([T^ . For 
p ^ the leading term of a considered equation of state was shown to be given by kpc^ = s (/cp)^°, where s is 
a dimensionless constant, cf. (|Tn|) . 

The 2-parameter family of orbits that converges to B2 can be conveniently characterized by the constants A 
and B, defined according to 

A := lim U{1 - 1/)-^° , B := lim l]/o(-o+i)(l - V)-i . (37) 
Expressed in these quantities we obtain the following, when r ^ 

p{r) = k-^B""^^^ (5r^°+i + i^-i^^o+i (^2(no + 1) + k'^ B''^^^ A'^ {'^ + 2no)^ (5r^°+' + O(-) 

(38a) 

m(r) = M(l - A(5r^°+i) + O(-) (38b) 
p{r) = ^MASr^^^O{-), (38c) 

where Sr = 

The radius R and mass M for a solution are uniquely determined by the values of A and B through 

" 4^ s^2B ' " 4^7g3 {s + 25)3 • (^9) 

The dimensionless quotient GM/{c^R) is given by GM/{c^R) = B/{s + 25), so that the line element (QJ on 
the surface r = is determined by e^^ = c'^e-^^ = 1/(1 - ^f^) = 1 + 2B/s. (Hence, naturally, GM/{c^R) 
can be used instead of B together with A to parameterize the different solutions.) 

Eq. ((33 and the related formulas are useful in, e.g., numerical computations, since it is fairly easy to compute 
A and B. Below, in numerical applications, we are going to compare solutions with a reference solution, i.e., we 
are going to consider the dimensionless ratios R/R^q^, M/M^q^^ where i?ref and Mref are the radius and mass 
of a typical reference solution. Of course, the dimensional factors (A:c^)/(47rsG) and {kc^){4:7rsG^) drop out in 
this case. 

For illustrative purposes consider cylindrical coordinates at B2, i.e., the coordinate transformation ([/, V, O) = 
(5sin^, 1 — ecos^, h^) with < ^ < 7r/2, < /i, < e, which leads to 

A= lim e^-^°sin^(cos^)-^° , 5 = lim (ecos^)"^ . (40) 

As every orbit in a neighborhood of B2 converges to ^2 as A ^ 00, we obtain approximate expressions with 
arbitrary accuracy by choosing e sufficiently small. In this picture, every orbit is uniquely characterized by its 
intersection point (/i, (j)) with the small e-cylinder, and accordingly these values uniquely determine M and R 
by combining ((5^ with (00)). A similar discussion holds for L4 when no = 0; for the corresponding Newtonian 
discussion, see jllj . 
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6 Mass-radius theorems 

In this section we will prove several theorems concerning mass-radius properties of solutions, where we focus on 
regular solutions. The underlying equations of state are as always understood to be asymptotically polytropic 
for p ^ and asymptotically linear for p ^ oo, as discussed in Sec.|31 

Theorem 6.1. All regular solutions have infinite masses and infinite radii if < | ctTid a < 1. 

Proof. Consider the function ^2(^5 V, 1^), defined by 

^2 := U{S -9V)^7V-6 . (41) 

The surface ^2=0 coincides with the regular orbit from L2 to Bi for no = 5 when projected onto the = 
plane, cf. ([T^ . Taking the derivative of this function and evaluating it on the surface ^2 = yields 

^ = V{1-V) {{5- 6T)(1 + 2 a ((4/7 + /7a - 1)(3 - 4/7) - 3T(1 - /7)(1 + Ua)) ) . (42) 

When Ftv < I (T > |) and a < 1, then d^2/d\ < 0. Since in addition the unstable subspaces of the fixed 
points L2 lie in the region ^2 < of the state space, it follows that a regular orbit can never leave this region. 
The only attractor in this part of the state space is the fixed point ^4, cf. Theorem 14. II Since B/^ gives rise to 
perfect fiuid solutions with infinite masses and radii the claim of the theorem is established. ■ 

Remark. The assumption a < 1 is just the dominant energy condition, and is satisfied, e.g., for causal relativistic 
polytropes. In the theorem, the condition a < 1 is a sufficient but not necessary; the statement is valid for 
much larger a, but not for arbitrarily large values. 

Let us define ai := (^-4 + v'16 + 2 (STat - 6) ^ /6. 

Theorem 6.2. All regular solutions with Vtc ^ (0, ^^max) have finite masses and radii ifT^ > | ctnd if a < ai 

on (0,l]max). 

Proof. Consider again the function ^2 as defined in If Fat > 6/5 (i.e., T < 5/6), then d^2/d\ > is 

satisfied, provided that a < ai. Under the same conditions, the unstable subspace of a fixed point on L2 lies 
in the region ^2 > 0- It follows that a regular orbit, which originates from a fixed point (3/4,0,r^c) on ^2, is 
confined to the region ^2 > of the state space, if a{Q) < ai for all Q < Qc- The only attractor in the region 
^2 > is the fixed point B2, and since B2 generates solutions with finite masses and radii, the theorem is 
established. ■ 

Remark. In it is shown that if 6/5 < Fq < 2 then the solution possesses a finite radius provided that a 
is sufficiently small (Theorem 4, p. 994). In [2H1 and ^01 theorems related to Theorem 16.21 are proved with 
completely different methods. In the above theorem, the lower limit for a is not a particularly good lower 
bound for "most" equations of state. Indeed, as follows from the next theorem, when no < 3, then all solutions 
have ffnite masses and radii for all Qc ^ [0, 1] and thus all values of a are allowed. 

Theorem 6.3. (Finiteness of perfect fluid solutions). All regular and non-regular perfect fluid solutions have 
finite radii and masses if uq < 3 (i.e., Fq > 3/4^. 

Proof. This theorem follows from Theorem 14.21 which shows that all solutions end at B2 if no < 3 (at I/4 if 
no = 0), and from that B2 {L/^) is associated with solutions with finite masses and radii. ■ 

Remark. {General validity). The last theorem holds irrespective of the asymptotic behavior of the equation of 
state at high pressures. Also for the first two theorems the asymptotic high pressure regime is irrelevant if one 
restricts the attention to solutions with finite central pressures. This is because the proofs in such cases do not 
rely on the inclusion of the boundary = 1 in the state space. This is in contrast to the next theorem which 
makes use of the Vt = 1 subset and relies on our assumptions about differentiability of the equation of state 
when Vt ^ 1 {p ^ 00). 
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Theorem 6.4. (Spiral structure of the {M, R)- diagram). Let < no < 3.^^ For sufficiently high central 
pressures the mass-radius diagram exhibits a spiral structure, i.e., {R{pc), M{pc)) is given by 

R{pc) \ _ ( Ro 



where Ro and Mo are constants, B is a matrix with positive determinant, and b a non-zero vector. The matrix 
J{^) G S0(2) describes a (positive) rotation by an angle (f>, and the constants ri and T2 are given by 

ri = — , r2 = — Vc , 44 

471 471 

where ^fc = + 4471 - j-^ see Table[^ 

Sketch of proof. The basic observation is that the regular orbit on 1] = 1 forms a spiral as it converges to T3, 
and that this spiral subsequently reflects itself as a spiral in the (i?, M)-diagram. 

Choose e small and set a;^ = 1 — e. The dynamical system asymptotically decouples as ^7 ^ 1 (recall from (|T^ 
and that we can always choose co in such a way that ^ = 0((1 — ^)^) and ^ = 0((1 — ^)^) as ^ 1), 
which enables us to approximately solve the system on [0,1]^ x [l]^,!]]: let C/l(A), ^(A) denote the regular 
solution of the system ((5a|) and ()8b|) with = l,cr = <Ji=7i — 1. Then (/7l(A), Vl{X), ^(A)) is an approximate 
solution of © when ^ = -/i(l - n)Hi^, where i^L = H{Ui^, Fl, !)• 

For large A the orbit (/7l(A), Vl{X)) has the form of a spiral, i.e., 

( vl[\) ) = ( ) +^M-Si\)B'J{52X)b', (45) 

where Ut^^ Vt^ are the coordinates of the fixed point T3; B' is a matrix with det > 0, and b' a non-zero vector; 
^1 = 4 (-^1^2) (3"^i ~ ^) ^2 = 4(^?^2) V^* "^^^ constants Si and 62 are the real and imaginary part of the 
complex eigenvalue that is associated with T3, cf. with Table Q 

Inserting ((l5j) into the equation for Q we obtain the approximation 

1 - n{X) = (1 - Qc) const exp {fiHr, (A - A)) , (46) 

where Ht^ = HCUts^Vt^^I) = 7i/(2 + 7^), and A and const are independent of ftc (the value of ^^(A) as 
A —00). 

Consider the orbit that originates from T3 into the state space: it intersects the plane Q = in di point / 
with approximately the same ([/, T/) -coordinates as T3. A regular solution with sufficiently close to 1 can be 
described by and It intersects 1] = at A = Ac, where 

Ac 2 — -^r + const = logpc + const . (47) 

(2-0 i — iZc -^Ts 

For the components ([/(Ac), l^(Ac)) of this regular solution we have {U (Xc) ^ V (Xc)) = {Uts^Vts) + (SUo^dVo) 
with 

<5C/o(??c) 



= (^) ' J(T2l0gPc)5", 



(48) 



where B^' (det S'' > 0) and b'^ have absorbed the constants. From ((l8j) we see that the set of regular solutions 
intersects 1^ = in a spiral in a neighborhood Uj of / in this plane. 

The orbit that originates from T3 and passes through I G Uj eventually intersects the mass-radius cylinder 
(cf. (|4()|) ) at some point O G Uo ^ cylinder, where the cylindrical coordinates of O determine the mass Mo and 
the radius Ro of the solution (cf. ((39j) ). The fiow of the dynamical system induces a diffeomorphism Uj Uo, 
so that the spiral in /7/ is mapped to a distorted (positively oriented) spiral in Uo- Since the cylindrical 
coordinates (0, h) in the small neighborhood Uo are related to (i?, M) by a positively oriented linear map, the 
spiral ((l8j) . with different B and 6, appears also in the (i?, M) diagram. ■ 

^-•^More generally, we can consider equations of state with arbitrary no as long as the orbit that originates from T3 ends at B2 
(or L4 if no = 0). 
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Remark. Note that the theorem does not rely on the behavior of the equation of state in the intermediate non- 
asymptotic regimes. Hence the theorem describes a universal phenomenon and connects it with the self-similar 
solution that corresponds to T3. A similar theorem has been proved in 24 using quite different methods. 
However, an advantage with the present approach, apart from less restrictive assumptions, brevity of proof, and 
clarity, is that it visually shows the importance of key solutions like the self-similar one associated with T3. 



7 Examples 

In this section the dynamical systems approach to relativistic stellar models will be illustrated by several exam- 
ples. The presented results also demonstrate the usefulness of the new framework in numerical computations. 



7.1 Relativistic polytropes 

As discussed in Sec. |2l a relativistic polytrope with index T can be described by the dynamical system Q 
with <j = (r — 1)Q and T = (1 + (F — l)Q)/r. The numerical integration of this system is straight forward, 
when we use the local analysis near the fixed points. The radius R and the total mass M of a solution can 
be obtained by various methods, in particular one can use equations (|35d|) and (|l^5e|) . or one can include the 
auxiliary equations (|35a|) in the numerical integration of the dynamical system; Eq. provides a different 
approach. Mass-radius diagrams for the regular solutions associated with different values of T are shown in 
Fig.ISl 

The theorems of Sec. El guarantee that all regular solutions possess finite R and M for no < 3 (F > 4/3). 
For 5 > no > 3 (F < 4/3), regular solutions with small central pressures must be finite fluid bodies. The 
numerical investigation improves these results: all regular solutions associated with F > 1.2971 (no ^ 3.366) 
possess finite R and M. A similar analysis was performed in [221 for purely polytropic equations of state: for a 
polytropic index n < 3.339 only finite fluid bodies appear. For an extensive discussion of the methods involved 
in establishing these results we refer the reader to that paper. 




R R R 

comp comp comp 

(a)r = 4/3 (b)r = 5/3 (c)r = 2.0 



Figure 5: Mass-radius diagrams for the regular relativistic polytropes. We display the compactified quantities 
^comp = (M/Mref)/(1 + M/Mref) and i?comp = (i^/i^ref )/(! + ^/^ref ) , whcrc M^ei and R^ei are typical values. 
The arrows indicate the direction of increasing i-e., increasing central pressures. 

Letting F ^ 1 for relativistic polytropes implies that <j ^ 0, and hence the dynamical system (|H|) approaches 
the Newtonian system, cf. Appendix H Moreover, the equation for Q decouples from the U- and V-equations, 
and the flows on the 1] = and the Q = 1 subset coincide in the limit. It follows that the regular surface 
is a surface perpendicular to the Q = const plane; its projection to 1] = const is depicted in Fig. |6(c)| Note 
also that the limit F ^ 1 reveals a relationship between the Newtonian self-similar solution associated with B4 
(discussed, e.g., in 21 ) and the relativistic self-similar solution associated with T3. 
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0.2 0.4 0.6 0.8 1 0.2 0.4 0.6 0.8 1 0.2 0.4 0.6 0.8 1 



(a)r = i.i (b)r = i.oi (c)r = 1.001 

Figure 6: Superimposed regular orbits for the relativistic polytrope, with varying F = Fq = 71. The dashed line 
corresponds to the orbit in the 1] = plane and the solid line to the orbit in the 1] = 1 plane. 

Regular solutions associated with an equation of state characterized by an arbitrary adiabatic index T{p) close 
to 1 lie close to the regular "soft-limit-surface" .^^ Hence, deviations from the soft-limit-surface measure both 
deviations of the equation of state from the soft limit and relativistic effects. Stiffer models are located "further 
out" compared to the soft-limit-surface and increasing stiffness implies that the [7 = 3/4 surface is approached, 
and it is on this "stiff-limit-surface" the incompressible perfect fluid solutions reside, as discussed previously. 
Note also that the Buchdahl surface, depicted in Fig. IHl for increasingly soft equations of state moves toward 
the V = 1 surface, thereby becoming increasingly less restrictive. This illustrates that the Buchdahl inequality 
is a purely relativistic effect, since the surface disappears from the interior state space in the soft, and thereby 
Newtonian, limit. 

7.2 The ideal neutron gas 

The equation of state of the ideal neutron gas provides a simple model for neutron star matter; it is given in 
implicit form in equation (|18|) . As discussed in Sec. [2l the ideal neutron gas can be naturally described in the 
dynamical systems framework, and is easy to handle numerically. 

In Fig. [3 regular orbits for various values of Qc and the regular orbits on ^7 = and Q = 1 are depicted. 

7.3 Composite equations of state 

Consider first a stiff perfect fluid characterized by the equation of state p = {p — p-)c^^ which can be rewritten 
as 

kpc^ = l^kp, (49) 

where k = (p_c^)~^. If we choose u = kp, we obtain / = 1, T = ft, and a = ft^ i.e., the dynamical system ^ 
becomes very simple. In this context, it is of interest to note that the solution that corresponds to the orbit 
from T3 to L4 (which is a solution with finite radius and mass but infinite central pressure) is explicitly known: 
it is a special case of a Tolman solution ^HI- the present formalism, the corresponding orbit is given by 
= 3/(5 + 1]) and 1/ = (2 + l])/(2 + 311 + 417^). (For interesting features of this solution, see 17 .) 

As an example of a composite equation of state, let us continuously join a stiff equation of state p = (p — p-) 
(when p > pj) to di relativistic polytrope (when p < pj)?^ Such an equation of state possesses a kink at pj. 

Consider an arbitrary one-parameter class of equations of state characterized by a parameter > 0, where Fq and 71 converge 
to 1 as /i ^ 0. Then the associated regular surface approaches the soft-limit-surface as ^ 0. 

^^If an equation of state is only known for low pressure (as is true in reality), then mass estimates for stellar models can be 
obtained by extending the equation of state as a stiff fluid, see e.g., 1 . 
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Figure 7: Typical regular orbits for the ideal neutron gas together with the orbits on the boundaries ^7 = and 



As outlined in Sec. [31 we work with two state spaces and dynamical systems, one pertaining to the relativistic 
polytrope and one representing the stiff fluid. The jump at pj reflects itself in a map between the two state 
spaces. Since p is continuous, the values of [/, V are continuous under the map, however, there is a jump in 
Since Preipoi(Pj) = Pstifffe) we have 



{kpj) 



i/r 



Pjj = p-c ^Pj . (50) 
Dividing by pj and using that uJreipoi{p) = {kp)^^~^'^^^ and a;stiff(p) = {P-c'^)~^Pj we obtain that 

^relpol = (^1) 

at the jump. If F = 2, i.e., if the relativistic polytrope is asymptotically stiff, then there is no jump at all. A 
typical solution with a jump is illustrated in Fig. [HI 

Note that, in general, if one starts with a regular solution in the first state space, then this solution has to 
be matched with a non-regular one associated with the other (extended) equation of state in the second state 
space. 



8 Concluding remarks 

In this paper we have derived a dynamical systems formulation for the study of spherically symmetric relativistic 
stellar models. The method of "homology invariants" known from the theory of Newtonian polytropes has been 
generalized both to a general relativistic context and to the broad class of asymptotically polytropic/linear 
equations of state. For previous related work see ^T] (Newtonian asymptotic polytropes) and [211 (polytropes 
in GR). 

The present formulation has turned out to be advantageous in many respects: (i) it provides a clear visual 
representation of the solution spaces associated with broad classes of equations of state; in particular it is 
revealed how the global qualitative properties of the solution space are influenced by the equation of state; (ii) 
the formulation makes the theory of dynamical systems available, which makes it possible to prove a number 
of theorems and describe the qualitative behavior of solutions; and (iii) the framework is particularly suited for 
numerical computations, since the numerics is supported by the local dynamical systems analysis. 

The idea of using dimensionless variables and exploiting asymptotic symmetries and properties to derive a 
dynamical systems formulation is a quite versatile one. As another example, in a future paper, we will show 
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Figure 8: An example of an orbit for a composite equation of state consisting of a relativistic poly trope (F = 
11/6) and a stiff fluid. We have chosen to cut the "stiff fluid orbit" at = 0.50 (where {U,V) = (0.58,0.60)). 
The jump transformation law then yields a map to the point ([/, V, = (0.58, 0.60, 0.75) in the state space of 
the relativistic poly trope, which serves as initial data for the subsequent evolution in that state space. 

how one also can treat a collisionless gas. The "dynamical systems approach" can even be applied to problems, 
in general relativity and other areas, when no symmetries exist at all as will be shown in another set of 
papers. Thus the ideas in the present work should be seen in a quite broad context, and it is likely that there 
are many problems in very different areas that could benefit from the type of ideas and the approach presented 
in this paper. 
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A Appendix: Comparison with the Newtonian theory 

The Newtonian equations of hydrostatic equilibrium are obtained from Q and Q by letting c ^ oo, however, 
note that in the Newtonian case p is interpreted as the rest-mass density and c^^ the Newtonian 

potential. Consistently, since a = p/{pc^)^ in the dynamical system we must set a = 0. Note also that v 
reduces io v = Gpm/{rp). Since a ^ for p ^ 0, the relativistic system (|H| and its Newtonian counterpart 
coincide on 1] = 0. 

In the Newtonian case we can treat equations of state p{p) with less restrictive asymptotic behavior in the 
regime p ^ oo. This is because it was the existence of a{Q) on the right hand side of the dynamical system that 
forced the assumption of asymptotic linearity upon us in the relativistic case. In particular we can naturally 
include asymptotically polytropic behavior as p ^ oo in the Newtonian case (cf. 11 ). 

In Newtonian gravity there exist two independent dimensional scales, space and time, while in relativity, space 
and time are connected by the speed of light c, so that only one single dimensional scale remains. A Newtonian 
equation of state can be written in the implicit form 



kip (/)(cj) , k2p x(^) , 



(52) 
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where c<; is a dimensionless variable and ki and k2 are constants carrying the dimension and respec- 
tively; and X monotone functions in cj containing any number of dimensionless parameters. An explicit 
representation of is k2p = ijjikip). 

Since our dynamical system is expressed in terms of purely dimensionless quantities, the dimensional parameters 
have to drop out. Consequently, a single dynamical system must be capable of describing the entire two- 
parameter set of equations of state ((5^ parameterized by ki and A:2. Indeed, the relevant function entering the 
dynamical system in the Newtonian case is T(a;) given by T{uo) = {d\nx/ d\iiix}){dhi(l) / dhiuo)~^ . Evidently, 
this function does not depend on ki and /c2. The same clearly holds for f{uo) = d\n<p/d\nu;. 

Example. As an illustrative example, consider the equation of state p = Cip^° + C2P^^ . By rescaling the 
equation and the two parameters one can write the equation of state as k2p = {kip)^° + (/cip)^^, and thus 
the equation has been brought to the desired form. The choice co = {kipY ensures that T(a;) and therefore 
also T(yt) are independent of ki and /c2, and interpolate monotonically between Tq and Ti. In analogy to 
the detailed discussion in Sec. [31 the constant a must be chosen sufficiently small so that a C^-differentiable 
dynamical system is obtained. 

Remark. {Consistent translatory invariance). As seen previously, it is the freedom to choose u;{p) that allows 
one to treat a one-parameter class of equations of state simultaneously by one specified dynamical system in the 
relativistic case. The freedom we were able to exploit in the Newtonian case to cover even two-parameter classes 
of equations of state reflects itself in the consistent translatory invariance of the Newtonian dynamical system. 
Namely, if x(^) (x(^) := {u{^)^ v{^)^u;{^)) (recall that ^ = Inr) is a solution of the Newtonian dynamical system, 
then so is xq := x(^ — ^o), and moreover, the translated solution gives rise to a perfect fluid solution associated 
with an equation of state with the same T{p). In contrast, in the relativistic case this consistency is broken: if 
x(^) is a solution of the dynamical system Q, then so is xo(^) = x(^ — ^o)- However, whereas the solution x(^) 
can be consistently interpreted as a relativistic perfect fluid solution, associated with an equation of state p{p) 
(determined by T and a), this is not the case with the translated solution. Although the translated solution 
satisfies differential equations associated with p(p), the initial data are not consistent with this equation of 
state. Accordingly, only one particular solution on every orbit of ©, can be interpreted as a relativistic perfect 
fluid solution. However, in the compactified dynamical system ((H| this 'defect' is remedied by the appearance 
of the freedom in the new independent variable A (this 'cure' could also have been implemented in the u^v^uo- 
formulation by defining ^ through dr jd^ = r, i.e, by letting ^ = Inr + ^o, where is an arbitrary constant, 
instead of setting ^ = Inr, however, the direct relation between ^ and r is sometimes useful). The reason behind 
the difference in the Newtonian and relativistic cases is the appearance of a in the relativistic dynamical system; 
cr uniquely determines the equation of state while T only does so up to a proportionality constant. 

To conclude the remarks on the Newtonian case, we consider the Newtonian counterpart of the mass-radius 
formulas ((S^ . These equations simplify considerably in the Newtonian case. Starting from the two-parameter 
family (|52|) we obtain for the surface potential ^^{K) = —GM/R the expression GM/R = k^^k2B/s, and for 
R^ and 
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